{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "# Design of vector-Apodizing Phase Plate patterns\n",
    "\n",
    "This tutorial introduces the basics of designing a phase pattern for a vector-Apodizing Phase Plate (vAPP) coronagraph.\n",
    "\n",
    "We'll start by importing all relevant libraries and setting up our pupil and focal grids. "
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "from hcipy import *\n",
    "import numpy as np\n",
    "import matplotlib.pyplot as plt\n",
    "from mpl_toolkits.axes_grid1 import make_axes_locatable\n",
    "\n",
    "pupil_grid = make_pupil_grid(512)\n",
    "focal_grid = make_focal_grid(4, 20)\n",
    "\n",
    "prop = FraunhoferPropagator(pupil_grid, focal_grid)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The design of the vAPP, i.e. the phase pattern, is uniquely defined by the pupil of a telescope or instrument, and the dark zone shape. Here we design a semi-realistic vAPP for the Magellan telescope using the existing make_magellan_aperture function."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "aperture = make_magellan_aperture(True)\n",
    "\n",
    "telescope_pupil = aperture(pupil_grid)\n",
    "\n",
    "imshow_field(telescope_pupil, cmap='gray')\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we are going to define the dark zone. For this tutorial, we would like the dark zone to have a 'D' shape and a contrast of $10^{-5}$ relative to the stellar core. \n",
    "\n",
    "We define the inner working angle by the separation from the star where the desired contrast level is reached and set it to of 2 $\\lambda/D$. The outer working angle is 15 $\\lambda/D$. \n",
    "\n",
    "The dark zone is created as a pixel mask on the focal grid. From this pixel mask we create a contrast map by multiplying with the contrast level."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "contrast_level = 1e-5\n",
    "\n",
    "dark_zone = (circular_aperture(30)(focal_grid)).astype(bool)*(focal_grid.x>2)\n",
    "\n",
    "contrast = focal_grid.ones()\n",
    "\n",
    "contrast[dark_zone] = contrast_level\n",
    "\n",
    "imshow_field(np.log10(contrast))\n",
    "plt.colorbar()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "We define a plot function for the vAPP that shows the phase pattern and PSF of a vAPP."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "def plot_vapp(vAPP, prop):\n",
    "    '''Plot the phase pattern and PSF of a vAPP\n",
    "\n",
    "    Parameters\n",
    "    ---------\n",
    "    vAPP : Wavefront\n",
    "        The wavefront of a vAPP mask, containing the vAPP pattern as phase and \n",
    "        the telescope pupil as amplitude\n",
    "    prop : Function\n",
    "        A propagator function that propagates the wavefront to a focal plane\n",
    "    '''    \n",
    "\n",
    "    # Plotting the phase pattern and the PSF\n",
    "    fig = plt.figure()\n",
    "    ax1 = fig.add_subplot(121)\n",
    "    im1 = imshow_field(vAPP.phase, mask=vAPP.amplitude, cmap='RdBu')\n",
    "\n",
    "    divider = make_axes_locatable(ax1)\n",
    "    cax = divider.append_axes('right', size='5%', pad=0.05)\n",
    "    fig.colorbar(im1, cax=cax, orientation='vertical')\n",
    "\n",
    "    ax2 = fig.add_subplot(122)\n",
    "    im2 = imshow_field(np.log10(prop(vAPP).intensity/np.max(prop(vAPP).intensity)),vmin = -5, cmap='inferno')\n",
    "\n",
    "    divider = make_axes_locatable(ax2)\n",
    "    cax = divider.append_axes('right', size='5%', pad=0.05)\n",
    "    fig.colorbar(im2, cax=cax, orientation='vertical')\n",
    "    plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now that everything is set up, we can start generating the vAPP phase pattern.\n",
    "\n",
    "We use an updated version of the Gerchberg-Saxton algorithm to do the calculation. This is an iterative method, so we set the number of iterations to be 80. \n",
    "\n",
    "We generate the input wavefront from the electric field in the pupil. Note that it is possible to use a pre-calculated wavefront that is optimized using different methods. \n",
    "\n",
    "The output wavefront called 'vAPP' has unity amplitude in the pupil and the desired phase pattern of the vAPP."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Setting up the vAPP calculation parameters.\n",
    "num_iterations = 80\n",
    "wavefront = Wavefront(telescope_pupil, 1)\n",
    "\n",
    "# Generate the vAPP pattern.\n",
    "vAPP = generate_app_keller(wavefront, prop, contrast, num_iterations, beta = 1)\n",
    "\n",
    "plot_vapp(vAPP, prop)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "For a classical apodizing phase plate, the design would be finished. \n",
    "\n",
    "However, the vector-apodizing phase plate applies gemeometric phase with a sign that depends on the circular polarization state of the incoming light. This is explained in detail in the \"VectorApodizingPhasePlate\" tutorial.\n",
    "\n",
    "In short, this means that without splitting the PSFs in circular polarization, the bright-side of one polarization state would be imaged on the dark zone of the orthogonal polarization state.\n",
    "\n",
    "This has three consequences. \n",
    "1. We need to separate the two polarization states by adding a phase ramp in the y-direction. This is the grating-vAPP. \n",
    "2. We need to add dark zones on different locations in the focal plane to make sure that the two coronagraphic PSFs do not add light in the dark zone of the orthogonal polarization state. \n",
    "3. A polarization leakage PSF, that is not apodized by the vAPP, will be imaged on the optical axis. This PSF can be used as a photometric reference if there is a dark zone for both coronagraphic PSFs there as well.\n",
    "\n",
    "We decide that the central circular dark zone has a radius of 5 $\\lambda/D$. Combined with the 15 $\\lambda/D$ radius of the previous dark zone, we decide that we need a separation of 20 $\\lambda/D$. \n",
    "\n",
    "Because we are adding the phase ramp to split the polarization states in a later stage, we calculate the locations of all three dark zones from the coordinate system of one coronagraphic PSF. \n",
    "\n",
    "Therefore, we do not need to shift the D-shaped dark zone to 20 $\\lambda/D$ in the y-direction, but leave that on-axis. We shift the circular dark zone by -20 $\\lambda/D$ in the y-direction, and add a second D-shaped dark zone with flipped orientation at -40 $\\lambda/D$ in the y-direction. A larger focal grid is needed to contain the new dark zone mask."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "# Updating focal grid\n",
    "focal_grid_upd = make_focal_grid(4, 60)\n",
    "\n",
    "# Generating the three dark zones\n",
    "dark_zone_upd = (circular_aperture(30)(focal_grid_upd)).astype(bool)*(focal_grid_upd.x>2)\n",
    "dark_zone_upd += circular_aperture(10)(focal_grid_upd.shifted((0,20))).astype(bool)\n",
    "dark_zone_upd +=  (circular_aperture(30)(focal_grid_upd.shifted((0,40)))).astype(bool)*(focal_grid_upd.x<-2)\n",
    "\n",
    "# Making the contrast map\n",
    "contrast = focal_grid_upd.ones()\n",
    "contrast[dark_zone_upd] = contrast_level\n",
    "\n",
    "# Plotting the contrast map\n",
    "imshow_field(np.log10(contrast))\n",
    "plt.colorbar()\n",
    "plt.show()"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "With this updated dark zone, we can calculate the phase pattern again using the same method as before.\n"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "prop2 = FraunhoferPropagator(pupil_grid, focal_grid_upd)\n",
    "\n",
    "num_iterations = 80\n",
    "contrast_level = 1e-5\n",
    "\n",
    "vAPP = generate_app_keller(wavefront, prop2, contrast, num_iterations, beta = 1)\n",
    "\n",
    "plot_vapp(vAPP, prop2)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now we add the phase ramp in the y-direction and the vAPP pattern is finished."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {},
   "outputs": [],
   "source": [
    "gvAPP = vAPP.electric_field * np.exp(1j * pupil_grid.y * 2 * np.pi * 20) * telescope_pupil\n",
    "\n",
    "gvAPP = Wavefront(gvAPP)\n",
    "\n",
    "plot_vapp(gvAPP, prop2)"
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Python 3",
   "language": "python",
   "name": "python3"
  },
  "language_info": {
   "codemirror_mode": {
    "name": "ipython",
    "version": 3
   },
   "file_extension": ".py",
   "mimetype": "text/x-python",
   "name": "python",
   "nbconvert_exporter": "python",
   "pygments_lexer": "ipython3",
   "version": "3.6.5"
  },
  "level": "intermediate",
  "thumbnail_figure_index": 2
 },
 "nbformat": 4,
 "nbformat_minor": 2
}
